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We construct a hexagonal lattice of repressing genes, such that each node represses three of 
the neighbors, and use it as a model for genetic regulation in spatially extended systems. Using 
symmetry arguments and stability analysis we argue that the repressor-lattice can be in a non- 
frustrated oscillating state with only three distinct phases. If the system size is not commensurate 
with three, oscillating solutions of several different phases are possible. As the strength of the 
interactions between the nodes increases, the system undergoes many transitions, breaking several 
symmetries. Eventually dynamical frustrated states appear, where the temporal evolution is chaotic, 
even though there are no built-in frustrations. Applications of the repressor-lattice to real biological 
systems are discussed. 
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Our understanding of genetic regulation inside the cell 
has greatly improved in recent years. A number of genetic 
circuits have been quantitatively characterized, ranging 
from switches to oscillators made up of negative feedback 
loops. The latter class of circuits is ubiquitous in regu- 
latory networks with oscillating gene expressions, two of 
the most important examples being the NFkB network 
for inflammatory response P, 0, Q and the p53-mdm2 
system which regulates cell apoptosis 0, Q . 

However, decisions taken inside the cell may depend 
crucially on the environment and may be cooperative, 
i.e. depend on the behavior of neighboring cells. This 
calls for theoretical modeling which explicitly takes the 
spatial arrangement of different cells into account. As a 
basic unit, we consider a negative feedback loop consist- 
ing of three proteins that repress each other by blocking 
the associated genes, which Leibler and Elowitz termed 
the 'repressilator' [6|. Previously, others have studied 
coupled repressilators to investigate quorum sensing Q 
and cell-to-cell communication [8|. As a further step, 
one might consider systems made up of regular arrays 
of cells interacting in a specific manner with neighbor- 
ing cells. Because of close packing, real arrays of cells in 
planar tissues often display hexagonal or near-h exag onal 
structure, e.g. in hepatic or retinal tissue [^, [13, 

Here we approach this general problem by extending 
the simple repressilator to a repressor-lattice - a hexag- 
onal array of repeated and overlapping repressilator mo- 
tifs, as shown in Fig. [1] Each node is repressed by three 
neighboring nodes and at the same time represses three 
other neighbors. A biological implementation of such a 
system would require a tissue where cells communicate 
specifically with their immediate neighbours, rather than 
in a mean-field manner as in quorum sensing. Such direct 
communication is in fact quite common, either through 
small conduits that connect the cells, or via proteins 
that span the membrane of the cells [I4I. Further, the 
directed nature of the interactions would require some 



form of epigenetic gene silencing, resulting in adjacent 
cells expressing different gene s even though they have 
exactly the same DNA [H, 14, 1^. The modeling frame- 



work we propose is, however, general and can be used to 
describe other kind of interations, such as bidirectional 
ones, which might be easier to realize experimentally. 

The lattice in Fig. [1] can be naturally constructed to 
be translationally invariant and such that all local loops 
are repressilator motifs. We will approach the problem 
by imposing periodic boundary conditions in order to 
preserve translational invariance. Later, we will discuss 
how these results translate to the case of a large lattice 
without the periodic boundary conditions, which is more 
relevant for biology. 
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FIG. 1: The construction of the repressor-lattice from 'units' 
of single repressilators suitably placed on a hexagonal lattice. 
Each link symbolizes a repressor between two nodes corre- 
sponding to repressing genes, proteins, species, etc. 



The basic repressilator motif may exhibit an oscillating 
state with a phase difference between consecutive vari- 
ables equal to 2tt/3. One can ask whether the entire 
repressor-lattice might exist in an oscillatory state where 
only three different phases are allowed, each differing by 
2tt/3. We will show that this is indeed the case, but 
lattice commensurability effects may break this scenario. 

In the repressor-lattice, the variable at a node (to, n) 
is repressed by three neighboring nodes which we repre- 
sented by an interaction term Fint , leading to a dynami- 
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FIG. 2: Systems of 3 x 3 (a), 4 x 4 (b) and 5 x 5 (c) nodes 
subjected to periodic boundary conditions as indicated by the 
extra links. The numbers refer to the different phases of the 
solutions just above Hopf bifurcations. In (a) the solution 
exhibits symmetry with respect to rotations of angles which 
are multiples of 2n/3. In (b), (c) this symmetry is broken, so 
that 3 distinct solutions coexist above the Hopf bifurcation. 
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cal equation for the concentration of species Xr, 



dt 



7 Xr. 
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We consider two types of interaction terms - either an 
additive repression (an 'or gate'): 
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or a multiplicative repression (an 'and gate'): 
111 
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In either case we use standard Michaelis-Menten terms 
to model the repression. The parameter c measures the 
constitutive production of the proteins, 7 determines the 
degradation rate and a the strength of the repression by 
another protein. Further, K is the dissociation constant 
of the binding complex whereas h is the Hill coefficient 
measuring its cooperativity. For simplicity we assign the 
same parameter values to all the nodes in the lattice. We 



note that Ref. [a| also introduced the associated mRNA 
for each gene resulting in six coupled ordinary differen- 
tial equations. For simplicity we keep only the protein 
variables leading to three coupled equations - a single re- 
pressilator with this sirnplification can still be brought 
into an oscillating state [16[ . 




FIG. 3: Solutions of repressor-lattices of sizes a): 3 x 3, b); 
4 X 4, c): 5 X 5, d): 6 X 6 with multiplicative interactions, 
Eq. ® and parameters c = 0.1, 7 = 1.0, K = 1.0, h = 2. The 
value of Q is in each case chosen to be just above the Hopf 
bifurcation. Note that three, four and five different phases 
exist in a), b), c), respectively. In d) there are however only 
three different phases. 



For a single repressilator there exists a large regime of 
parameter space where oscillations are possible [y, [Tij . 
The transition to oscillations occurs via a Hopf bifurca- 
tion. We find similar behavior in the repressor-lattice. 
As a starting point, a lattice with 3x3 nodes as in Fig. 
[2^ was simulated both with additive, Eq. ([2]), and mul- 
tiplicative, Eq. ([3]), couplings. Just above the Hopf bi- 
furcations, we found smooth oscillations with only three 
distinct phases as indicated by the numbers 1,2,3 labeling 
the nodes in Fig. [5^. The oscillating time series is shown 
in Fig. [3^. These solutions are trivially related to the 
solutions of the basic repressilator motif since each node 
receives three identical inputs, with a 27r/3 phase shift 
with respect to itself. Note that this solution is invariant 
under lattice rotation of multiples of 27r/3 

However, this scenario is not completely general. For 
instance, in the case of a lattice of size 4x4 (16 nodes) 
the corresponding dynamical solutions are different, as 
shown in Fig. [SJa. As in the previous case (Fig. [3^) we 
are relatively close to the first Hopf bifurcation. However, 
now phases of the oscillating solutions differ by 27r/4 be- 
tween the nodes. The origin of this is a commensurabil- 
ity effect between the number of nodes in the lattice and 
the associated number of possible phases of the oscillat- 
ing solutions. This commensurability effect is of course 
enforced by the periodic boundary conditions. The com- 
plete structure of the phases is shown on the 4x4 unit 
cell in Fig. ^Bp. The case of 5 x 5 is also affected by 
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commensur ability effects, as shown in Figs. [2j; and[3j:. 

As opposed to the 3x3 system, here the inputs ar- 
riving to a specific node are different. This reflects the 
fact that the oscillatory solution is no longer rotationally 
invariant. We note that all lattices which are commen- 
surate by three, i.e. 6x6 (see Fig. [5Ji), 9x9, etc, 
allow a non-frustrated, symmetric state similar to the 3 
x 3 system. These periodic solutions all exhibit a Gold- 
stone mode in the sense that it is possible to slide the 
phases as long as the phase differences are kept constant. 
This means that the specific values of the phases for the 
solutions arc determined by the initial conditions. 

In order to understand these solutions in depth, we 
perform a stability analysis. We consider the "or" gate 
Eqs. I1|2I and, since the system is translationally invari- 
ant, we search for a constant homogeneous solution: 
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{-1X*-c){K''+x*'''). (4) 



The equation for x* always has one, and only one, real 
positive solution. The next step is to perturb the homo- 
geneous solution in order to perform a stability analysis. 
We consider a first order perturbation of the form: 
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Notice that since the solution must have the periodicity 
of the lattice, km and kn should be natural numbers and 
also 1 < km,kn < L. Plugging the solution into Eq. 
P]) and expanding to first order in e yields the following 
dispersion relation: 
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where a = ahK^ {x*)^"^ / [K^' + {x*f]^. Other kind of 
interaction terms lead to the same dispersion relation 
simply with a slightly different definition of x* and a; 
for example, taking multiplicative interactions leads to 
h = a{h/K){x* lKf-'^/[l + {x*/K)^Y'. Eigenvalues A 
with a positive real part will destabilize the homogeneous 
solution. Taking the real part of expression ^ , the eigen- 
value with the largest real part is the one that minimizes 
the function: 
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Before finding the solutions, we stress that f{km,kn) = 
fi—krn, —kn), whilc the imaginary part of the eigenvalue 
changes sign when the wave vector changes sign. This 
means that the two vectors {km,kn) and {—km,—kn) 
minimizing the function / are the complex conjugate pair 
that will cause the Hopf bifurcation. The function / is 
independent of the parameters of the system, meaning 
that the kind of pattern depends not on the form of the 



interaction (as long as the lattice is homogeneous with 
the same geometry), but on the number of sites in the 
lattice. The value of 7 determines only how much we 
have to increase a to encounter the Hopf bifurcation. 

Plotting the function f{km/L,kn/L) in the first 
Brillouin zone, < fcm,A:„ < L we see that it 
achieves its absolute minimum for the couple of eigen- 
values {km/L,kn/L) = (1/3,2/3) and (km/L,kn/L) = 
(2/3,1/3), where f(km,kn) = -3/2, see Fig. H This 
means that a Hopf bifurcation will occur when 3a = 27. 
Of course these wave vectors are allowed only when L is 
a multiple of 3, so that the values of km and fc„ at the 
minimum arc natural numbers. 
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FIG. 4: The landscape of the function Eq. ((Tjl. Red dots mark 
absolute minima, corresponding to the symmetric solution for 
L multiple of 3. Blue lines mark the bottom of the valleys 
around the minima. When L is not multiple of 3, the absolute 
minima are not achievable, and the degenerate solutions are 
given by 3 complex conjugate pairs along these valleys. 

We can of course minimize the function / also for val- 
ues of L that are not multiples of 3. For L = 4, the 
minimum is /(I, 3) = /(3, 1) = —1, but also /(2,3) = 
/(3,2) = -1 and /(3,4) /(4,3) = -1. The case 
L = 5 is also a degenerate case. The minimum is 
/(2, 3) = /(3, 2) « -1.30902, but also /(I, 3) = /(4, 2) « 
-1.30902 and /(2,4) = /(3, 1) « -1.30902. Ah cases 
that are not multiples of 3 have this degeneracy, due to 
the symmetry of the lattice. Close to the Hopf bifur- 
cation, the number of observed phases will reflect the 
periodicity of the eigenfunction. In particular, there will 
always be 3 distinct phases if L is a multiple of 3 and L 
phases if L is a prime number. 

The phases of the eigenfunctions can be used to flg- 
ure out how the oscillation pattern will look like on the 
lattice: sites on the lattice at a distance Am, An such 
that fcmAm + fc„An ~ will be in phase. Fig. [2t shows 
one of these solutions of the 5x5 lattice, namely /(4, 1). 
The other 'symmetric' solutions can be obtained through 
rotations of multiples of 27r/3, respecting the hexagonal 
rotational symmetry of the lattice. 

One might expect that the symmetric solutions with 
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FIG. 5: A specific chaotic solution for a 5 x 5 lattice with 
multiplicative coupling (Lyapunov exponent equal to 0.028) 
obtained at a coupling strength equal to a = 12.8 with the 
variable in node (4, 5) plotted against the variable in node 
(2,2). Other parameter values are: c — 0.1,7 — 1.0, K — 
1.0, h — 3. Inset: Bifurcation diagram, showing the maxima 
and minima for dynamical solutions in the node point (4, 5) 
after a transient period of 15000 time units. The extremal val- 
ues are plotted against varying values of the coupling strength 
a. The critical value for the Hopf bifurcation calculated from 
the dispersion relation |S] (see main text) gives Oc = 1.838. 



boundary conditions to represent an isolated planar tis- 
sue, in which cells at the boundary receive no external 
signal. Wc found from simulations that such a system 
shows frustration effects similar to the case of periodic 
boundaries: when the steady state is destabilized, cells 
far from the boundaries exhibit the three-phase dynamics 
of the repressilator circuit, while closer to the boundaries 
the dynamics is more irregular, with more phases possi- 
ble. We did not observe any chaos in this case, even for 
very large values of a. 

In conclusion, the lattice model we have investigated 
here provides a simple starting point to study regulation 
in spatially extended biological systems. Future direc- 
tion could include, for instance, introducing an intrinsic 
'frustration' in the repressor-lattice. There are several 
ways of doing this, e.g. by lattice defects, or mutations 
modifying some of the interactions. For example, one 
can consider what happens when a specific repressor link 
is mutated into an activator. These generalizations may 
provide a useful framework for describing more specific 
cases of cell-to-cell communication in biological tissues. 

We are grateful to Namiko Mitarai, Joachim Math- 
iesen and Szabolcs Semsey for discussions. This work is 
supported by Danish National Research Foundation. 



five different phases. Fig. could exist even when pa- 
rameters are varied. This is not the case: when the cou- 
pling parameter a is increased, several transitions related 
to strong non-linear effects take place. Starting out with 
smooth periodic solutions of five distinct phases, am- 
plitude modulations set in when the coupling constant 
a w 2.6. At the same parameter value, we also observe 
that some phases that were distinct before this transition 
now coalesce with each other. At higher a- values, am- 
plitude modulations become even more pronounced and 
furthermore temporal period-doubling sets in. Increas- 
ing a even more, chaotic solutions eventually appear as 
shown in the bifurcation diagram and the attractor of 
Fig. [5] (similar bifurcation diagrams are observed for re- 
pressive cell-cell communication Q). Even though the 
lattice is made up of simple rcprcssilators without local 
frustration, the resulting dynamics is chaotic: it is not 
possible to keep the simple five-phase solutions when the 
repressilators arc coupled strongly with the neighbors. 
Each node in the 5x5 lattice exhibits a different bifurca- 
tion diagram (no simple period five symmetry operations 
are present) showing that all symmetries are eventually 
broken through a series of non- linear transitions. 

One may wonder how realistic periodic boundary con- 
ditions are for modeling real biological systems. In the 
3x3 case, this might be implemented in a single cell 
with 9 different genes, each repressed by three different 
ones. Having in mind extended systems, a more realistic 
case is to consider a large, finite lattice with non-periodic 
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